

In [1]:

    
%pylab inline

import pylab as pl

import numpy as np 
import matplotlib.pyplot as plt

import math


rmax = 0.8333333 # maximum radius of the suboff
xb = 3.333333     # bow
xm = 10.645833    # parallel plate 
xa = 13.979167    # afterbody  
xc = 14.291667    # after cap
cb1 = 1.126395101 
cb2 = 0.442874707
cb3 = 1.0/2.1
rh = 0.1175
k0 = 10.0
k1 = 44.6244
iDx_offset = 0.01          # offset 
idx = 1000# number of points
idy = 1000                 # number of points

total_length = 14.291666

#***************************************************************************
#                         CREATE THE CURVATURE OF THE BOW
#***************************************************************************

x_le_bow = 0.0
x_tr_bow = 3.333333
#x_tr_bow = 14.291666


x_bow = np.zeros(idx + 1)
y_bow = np.zeros(idx + 1)
a_bow = np.zeros(idx + 1)
b_bow = np.zeros(idx + 1)
fun = np.zeros(idx+1)
#boundary conditions
x_bow[0] = x_le_bow
a_bow[0] = -1.0
b_bow[0] = 1.0
y_bow[0] = 0.0


for i in range(1,idx + 1):
    x_bow[i] = x_bow[i-1] + x_tr_bow/idx
    a_bow[i] = 0.3 * x_bow[i] - 1.0
    b_bow[i] = 1.2 * x_bow[i] + 1.0
    y_bow[i] = rmax*(cb1 * x_bow[i] * a_bow[i]**4.0 + cb2 * x_bow[i]**2.0 * a_bow[i]**3.0 + 1.0 - (a_bow[i]**4 * b_bow[i]))**(1.0/2.1)
        #bow[i,0] = x_bow[i]
        #bow[i,1] = y_bow[i]
   # if 3.333333<=x_bow[i]<=10.645833:
   #     y_parralel[i] = rmax
        #print y_parralel[i]
    #if 10.645833<=x[i]<=

#***************************************************************************
#                         CREATE THE PARALLEL BODY
#***************************************************************************  


x_le_parallel_midBody = 3.333333

x_tr_parallel_midBody = 10.645833-3.333333

idx_parallel = 50

x_parallel_midBody = np.zeros(idx_parallel+1)
y_parallel_midBody = np.zeros(idx_parallel+1)

# boundary conditions 

x_parallel_midBody[0] = 3.333333
y_parallel_midBody[0] = rmax 

for i in range(1,len(x_parallel_midBody)):
    x_parallel_midBody[i] = x_parallel_midBody[i-1] + x_tr_parallel_midBody/idx_parallel
    y_parallel_midBody[i] = rmax
    
# store the vector + the the offset 3.333333    

#***************************************************************************
#                         CREATE THE AFTER BODY 
#***************************************************************************
idx_after_body = 1000

x_le_after_body = 10.645833

x_tr_after_body = 13.979167 - x_le_after_body

x_after_body = np.zeros(idx_after_body+ 1)


y_after_body = np.zeros(idx_after_body + 1)



#boundary conditions 
x_after_body[0] = x_le_after_body
y_after_body[0] = rmax

# construct the axial coordinate
for i in range(1,len(x_after_body)):
    x_after_body[i] = x_after_body[i-1] + x_tr_after_body/idx_after_body

# constants
rh= 0.1175
k0 = 10.0
kl = 44.6244
a_0 = rh**2
a_1 = rh * k0

a_2 =  20.0      -     (20.0 * math.pow(rh,2.0))        -   (4.0 * rh * k0)    -   ( 0.333333 * kl)

a_3 = -45.0      +     (45.0 * math.pow(rh,2.0))        +   (6.0 * rh * k0)    +    (kl)

a_4 = ( 36.0      -     (36.0 * math.pow(rh,2.0))        -   (4.0 * rh * k0)    -    (kl)  )

a_5 = -10.0      +     (10.0 * math.pow(rh,2.0))        +   (rh * k0)          +    (0.333333 * kl)


# assign the xi function   and the first entry in the vector  
xi = np.zeros(idx_after_body + 1)    
#xi[0] = (13.979167-x_after_body[0])/3.333333
#calculate the y_afterbody functions 
#print x_after_body
for i in range(0,len(x_after_body)):
    xi[i] = (13.979167-x_after_body[i])/3.333333
    
for i in range(0,len(y_after_body)):
    y_after_body[i] = rmax*(a_0 + a_1*xi[i]**2.0 + a_2*xi[i]**3.0 + a_3*xi[i]**4.0 + a_4*xi[i]**5.0 + a_5*xi[i]**6.0 )**1.0/2.0

#***************************************************************************
#                         CREATE THE AFTERBODY CAP
#***************************************************************************
idx_after_body_cap = 100
x_le_afterbody_cap = 13.979167
x_tr_afterbody_cap = 14.291666 - x_le_afterbody_cap

x_after_body_cap  = np.zeros(idx_after_body_cap+1)
y_after_body_cap  = np.zeros(idx_after_body_cap+1)

# boundary conditions 

x_after_body_cap[0] = x_le_afterbody_cap

for i in range(1,len(x_after_body_cap)):
    x_after_body_cap[i] = x_after_body_cap[i-1] + x_tr_afterbody_cap/idx_after_body_cap

for i in range(0,len(y_after_body_cap)):
    y_after_body_cap[i] = 0.1175 * rmax * (1.0 - (3.2 * x_after_body_cap[i] - 44.733333)**2.0)**0.5


fig = plt.figure(figsize=(30, 2),dpi=1200, facecolor='w', edgecolor='k')
plt.plot(x_bow,y_bow,'x',linestyle='-',linewidth = '2',marker='o',markersize=5)
plt.plot(x_parallel_midBody,y_parallel_midBody,'x',linestyle='-',linewidth = '2',marker='o',markersize=5)
plt.plot(x_after_body,2.0*y_after_body,'x',linestyle='-',linewidth = '2',marker='o',markersize=5)
plt.plot(x_after_body_cap,2.0*y_after_body_cap/17.0212765957  ,'x',linestyle='-',linewidth = '2',marker='o',markersize=5)

plt.show()









    



Populating the interactive namespace from numpy and matplotlib

**PYTHON CODE TO GENERATE THE DARPA SUBOFF GEOMETRY

ASSYMETRIC HULL

The model is in ft. with a scale 1/24

Can be converted to meters

INCLUDED ARE:

BOW EQ. FOR 0.0 FT <- X <- 3.333333 FT,

PARALLEL MID-BODY EQ. FOR 3.333333 FT <- X <- 10.645833 FT,

AFTERBODY EQ. FOR 10.645833 FT <- X <- 13.979167 FT

AFTERBODY CAP EQ. FOR 13-979167 FT <- X <- 14.291667 FT.

Offsets are computed every 0.1 except in the first 0.5 ft. where they are at 0.01 ft

Initialize the enviroment

**PYTHON CODE TO GENERATE THE DARPA SUBOFF GEOMETRY

** ASSYMETRIC HULL **

The model is in ft. with a scale 1/24

Can be converted to meters

INCLUDED ARE:

BOW EQ. FOR 0.0 FT <- X <- 3.333333 FT,

PARALLEL MID-BODY EQ. FOR 3.333333 FT <- X <- 10.645833 FT,

AFTERBODY EQ. FOR 10.645833 FT <- X <- 13.979167 FT

AFTERBODY CAP EQ. FOR 13-979167 FT <- X <- 14.291667 FT.

Offsets are computed every 0.1 except in the first 0.5 ft. where they are at 0.01 ft

Initialize the enviroment

ASSYMETRIC HULL